package javavis.jip3d.functions;

import java.util.ArrayList;

import cern.colt.matrix.DoubleMatrix2D;
import cern.colt.matrix.impl.DenseDoubleMatrix2D;
import cern.colt.matrix.linalg.SingularValueDecomposition;

import javavis.base.Function3DGroup;
import javavis.base.JIPException;
import javavis.base.ParamType;
import javavis.jip3d.geom.Neuron3D;
import javavis.jip3d.geom.Normal3D;
import javavis.jip3d.geom.Plane3D;
import javavis.jip3d.geom.Point3D;
import javavis.jip3d.geom.Vector3D;
import javavis.jip3d.gui.Function3D;
import javavis.jip3d.gui.FunctionParam;
import javavis.jip3d.gui.ScreenData;
import javavis.jip3d.gui.ScreenOptions;
import javavis.jip3d.gui.dataobjects.PlaneSet3D;
import javavis.jip3d.gui.dataobjects.PointSet3D;

public class FPPatchSR4K extends Function3D {

	/**
	 * @uml.property  name="relacion"
	 */
	double relacion = 2/Math.sqrt(3);
	/**
	 * @uml.property  name="minRad"
	 */
	double minRad = 0.01;

	public FPPatchSR4K()
	{
		super();
		this.allowed_input = ScreenOptions.tNEURONSET3D;
		this.group = Function3DGroup.Model3D;
		
		FunctionParam p5 = new FunctionParam("Process Type", ParamType.LIST);
		String []paux = new String[2];
		paux[0] = "Lineal";
		paux[1] = "MultiThread";
		//paux[2] = "JCuda";
		p5.setValue(paux);
		FunctionParam p1 = new FunctionParam("Neighbor gen", ParamType.INT);
		p1.setValue(2);
		FunctionParam p2 = new FunctionParam("PlaneThick", ParamType.FLOAT);
		p2.setValue(0.05);
		FunctionParam p3 = new FunctionParam("Check integrity", ParamType.BOOL);
		p3.setValue(false);
		FunctionParam p4 = new FunctionParam("Source PointSet", ParamType.SCRDATA);
		
		this.addParam(p5);
		this.addParam(p1);
		this.addParam(p2);
		this.addParam(p3);
		this.addParam(p4);
	}
	@Override
	public void proccessData(ScreenData scrData) throws JIPException {
		result_list = new ArrayList<ScreenData>();
		PlaneSet3D new_plane_set = new PlaneSet3D(new ScreenOptions());
		new_plane_set.name = "patches" + scrData.name.substring(3);
		String processType = this.paramValueString("Process Type");
System.out.println("tipo: "+processType);
		int gen = this.paramValueInt("Neighbor gen");
		float thickness = (float)this.paramValueReal("PlaneThick");
		Object []elements;
		Neuron3D element;
		NormalFinder[] finders;
		Plane3D plane;
		
		PointSet3D pointSet = null;
		boolean checkIntegrity = this.paramValueBool("Check integrity");
		if(checkIntegrity)
			pointSet = (PointSet3D)this.paramValueScrData("Source PointSet");
		elements = scrData.elements();
		finders = new NormalFinder[elements.length];
		int cont;
		
		//compute normals
		for(cont=0;cont<elements.length;cont++)
		{
			element = (Neuron3D)elements[cont];
			if(checkIntegrity)
			{
				if(check(element, pointSet, 0.05))
					finders[cont] = new NormalFinder(element, thickness, gen, cont);
				
			}
			else
				finders[cont] = new NormalFinder(element, thickness, gen, cont);
		}
		
		//gather results
		try {
			for(cont=0;cont<elements.length;cont++)
			{
				if(finders[cont]!=null)
				{
					finders[cont].runner.join();
					if(finders[cont].result!=null)
					{
						plane = new Plane3D(new Normal3D((Point3D)elements[cont], new Vector3D(finders[cont].result[0],
								finders[cont].result[1], finders[cont].result[2])));
						plane.radius = finders[cont].result[3];
						new_plane_set.insert(plane);
					}
				}
			}
		} catch (InterruptedException e) {
			// TODO Auto-generated catch block
			e.printStackTrace();
		}
		result_list.add(new_plane_set);
	}
	
	private boolean check(Neuron3D neuron, PointSet3D points, double rad)
	{
		boolean ret = false;
		Object []neighbors = points.range(neuron, rad);
		if(neighbors.length>0) ret = true;
		return ret;
	}

	/**
	 * @author  miguel
	 */
	public class NormalFinder implements Runnable
	{
		Thread runner;
		/**
		 * @uml.property  name="origin"
		 * @uml.associationEnd  
		 */
		Neuron3D origin;
		int gen;
		double[] result;
		float thick;
		
		public NormalFinder(Neuron3D o, float t,int d, int n)
		{
			runner = new Thread(this, "thread"+n);
			origin = o;
			gen = d;
			thick = t;
			result = null;
			runner.start();
		}
		//@Override
		public void run()
		{
			int cont;
			int generation = 0;
			Neuron3D []vector;
			Neuron3D neurona;
			ArrayList<Neuron3D> neighbors = new ArrayList<Neuron3D>();
			ArrayList<Neuron3D> nextgen = new ArrayList<Neuron3D>();
			ArrayList<Neuron3D> visited = new ArrayList<Neuron3D>();
			ArrayList<Point3D> vectores = new ArrayList<Point3D>();

			double [][]matriz = new double [3][3];
			double []singular_values;
			DenseDoubleMatrix2D matrix;
			DoubleMatrix2D U;
			SingularValueDecomposition SVD;
			int vsize;
			double []datos;
			double minDist, dist;
			double sigmaN, sigmaT;
			double thickness;
			
			//neighbor searching
			vector = origin.getNeuronNeighbors();
			visited.add(origin);
			
			minDist = 1000;
			for(Neuron3D n: vector)	neighbors.add(n);
			do {
				while(!neighbors.isEmpty()) 
				{
					neurona = neighbors.remove(0);
					dist = origin.getDistance(neurona);
					if(dist<minDist) minDist = dist;
					vectores.add(neurona.subPoint(origin));
					visited.add(neurona);
					for(Neuron3D n: neurona.getNeuronNeighbors())
					{
						if(neighbors.indexOf(n)==-1 && visited.indexOf(n)==-1) nextgen.add(n);
					}
				}
				neighbors.addAll(nextgen);
				nextgen.clear();
				generation++;
			} while (generation < gen);
			
			//normal computing
			vsize = vectores.size();
			if(vsize>9)
			{
				matriz[0][0] = matriz[0][1] = matriz[0][2] =
					matriz[1][0] = matriz[1][1] = matriz[1][2] =
					matriz[2][0] = matriz[2][1] = matriz[2][2] = 0;
				for(cont=0;cont<vsize;cont++)
				{
					datos = vectores.get(cont).getCoords();
					matriz[0][0] += datos[0] * datos[0];
					matriz[0][1] += datos[0] * datos[1];
					matriz[0][2] += datos[0] * datos[2];
					matriz[1][1] += datos[1] * datos[1];
					matriz[1][2] += datos[1] * datos[2];
					matriz[2][2] += datos[2] * datos[2];
				}
				matriz[1][0] = matriz[0][1];
				matriz[2][0] = matriz[0][2];
				matriz[2][1] = matriz[1][2];

				matrix = new DenseDoubleMatrix2D(matriz);
				SVD = new SingularValueDecomposition(matrix);

				singular_values = SVD.getSingularValues();
				U = SVD.getU();
				singular_values = SVD.getSingularValues();

				sigmaN = singular_values[2];
				sigmaT = Math.sqrt(singular_values[0]*singular_values[1]);

				thickness = Math.atan(relacion * sigmaN/sigmaT);
				if(thickness<thick)
				{
					result = new double[4];
					result[0] = U.getQuick(0, 2);
					result[1] = U.getQuick(1, 2);
					result[2] = U.getQuick(2, 2);
					result[3] = minDist<minRad?minRad:minDist;
				}
			}
		}
	}

}
